************************************************************
************************************************************
*                                                          *
*      Part 3 : Empirical Analyses                         *
*                                                          *
************************************************************
************************************************************

cd ".\replication_folder"


****** 1. Summary graphs and statistics ********************


*** 1.1 Figure 1 Coal Production and Employment from 2001 to 2020

use merged_base, clear

collapse  (sum) n_close n_active production average_employees labor_hours , by( year)

tsset year

foreach v in  n_close n_active production average_employees labor_hours  {
	g `v'_2001 = `v' if year==2001
	carryforward `v'_2001, replace
	g `v'_index = 100*`v'/`v'_2001
}
 
tsline n_active_index production_index average_employees_index labor_hours_index if year>=2001, scheme(s2mono)  tline(2011,lcolor(gray))  xlabel(,grid)  ylabel( ,  angle(horizontal)) ylabel(, format(%9.0fc))  xtitle("") ytitle("")  ///
legend(on row(2) label(1 "Number of Mines") label( 2 "Coal Production") label( 3 "Average Employees") label( 4 "Labor Hours"))   title("")
graph export fig1.tif, replace


*** 1.2 Treatment and control group figures 

use merged_base, clear

g county=fips

* Figure 2 Treatment and Baseline Control Groups
maptile treat_oil_gas if year==2011 , geo(county2010) twopt(legend(rows(1) label(3 "Treatment") label(2 "Control") label(1 "Neither")) title("")) stateoutline(vthin)  savegraph(f2)  propcolor  replace
grc1leg  f2.gph ,scheme(s2mono)
graph export fig2.tif, replace

* Figure 3 Alternative Control Groups
maptile treat if year==2011 , geo(county2010) twopt(legend(rows(1) label(3 "Treatment") label(2 "Control") label(1 "Neither")) title("Panel A: Within Coal States")) stateoutline(vthin)   savegraph(f31) replace
maptile treat_100 if year==2011 , geo(county2010) twopt(legend(rows(1) label(3 "Treatment") label(2 "Control") label(1 "Neither")) title("Panel B: Within 100 Miles")) stateoutline(vthin)    propcolor savegraph(f32) replace
maptile treat_200 if year==2011 , geo(county2010) twopt(legend(rows(1) label(3 "Treatment") label(2 "Control") label(1 "Neither")) title("Panel C: Within 200 Miles")) stateoutline(vthin)    propcolor savegraph(f33) replace
maptile treat_psm if year==2011  , geo(county2010) twopt(legend(label(3 "Treatment") label(2 "Control") label(1 "Neither")) title("Panel D: PSM")) stateoutline(vthin)   propcolor savegraph(f34) replace
grc1leg  f31.gph f32.gph f33.gph f34.gph,scheme(s2mono)
graph export fig3.tif, replace


*** 1.3 Figure 4 Access to Finance and Mobility

use merged_base, clear

g county=fips
local s "treat_oil_gas"

maptile low_atf_bra if year==2011 & !missing(`s'), geo(county2010) twopt(legend(rows(1) label(3 "Yes") label(2 "No") label(1 "Neither")) title("Panel A: Low A2F-Banks")) stateoutline(vthin)  savegraph(f41) replace propcolor
maptile low_large_bra if year==2011 & !missing(`s'), geo(county2010) twopt( legend(  rows(1) label(3 "Yes") label(2 "No") label(1 "Neither")) title("Panel B: Low A2F-Large Banks")) stateoutline(vthin)  savegraph(f42) replace propcolor
maptile rural if year==2011 & !missing(`s'), geo(county2010) twopt(legend(rows(1) label(3 "Yes") label(2 "No") label(1 "Neither")) title("Panel C: Rural")) stateoutline(vthin)  savegraph(f43) replace propcolor
maptile low_mob_rel if year==2011 & !missing(`s'), geo(county2010) twopt(legend(rows(1) label(3 "Yes") label(2 "No") label(1 "Neither")) title("Panel D: Low Mobility")) stateoutline(vthin)  savegraph(f44) replace propcolor
grc1leg  f41.gph f42.gph f43.gph f44.gph,scheme(s2mono)
graph export fig4.tif, replace


*** 1.4  Summary statistics

use merged_base, clear

* Table A1 Correlations of Access to Finance and Economic Mobility
correlate  low_atf_bra low_large_bra rural low_mob_rel if year==2011 & !mi(treat_oil_gas)

* Table 1 County Characteristics of Treatment and Control
* Panel A: County Characteristics as of 2011
pstest emp_all avg_pay_all popestimate pop_old pop_young pop_white pop_black application_amt_1  low_atf_bra low_large_bra rural low_mob_rel  if year==2011 , raw t(treat_oil_gas)
* Panel B: Five-Year Average Characterize between 2007 and 2011 
pstest pop_tot black_per median_age median_inc median_home_value inequality single_parent education  if year==2011  , raw t(treat_oil_gas)



****** 2 Regressions ************************************************

global options "bdec(2) sdec(2) nor2 "


*** 2.1 Genertae Table 2, Table 7, and Table A3 

local outcomelist " emp_all avg_pay_all "                                                    // Baseline outcome variables for Tables 2, 7, and A3

*** 2.1.1 Table 2 Baseline Regressions: Private Sector Employment and Wages

local s "treat_oil_gas"                                                                      // Baseline treatment indicator for Tables 2

qui foreach outcome in `outcomelist' {

use merged_base, clear

replace share=. if share==0

egen i1 = xtile(share) if year==2011,  n(4)
g i2 = (i1==4) if !missing(i1)
egen share_dummy = mean(i2) , by(fips)

keep if tin(2007,2016)
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

g treat_post = `s'*post

* With only year and county FEs (e.g., column 1 of Table 2)
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table2.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

* Add interaction between each county characteristic and year FE (e.g., columns 2 to 5 of Table 2)
foreach c in median_inc median_home_value education black_per {
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._`c' ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._`c') 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table2.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC
}

* Account for all controls (e.g., column 6 of Table 2)
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._median_inc i.year#c._median_home_value i.year#c._education i.year#c._black_per) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._median_inc i.year#c._median_home_value i.year#c._education i.year#c._black_per) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table2.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}


*** 2.1.2 Table 7 Treatment Intensity Specification

local s  "share"                                                                           // Alternative intensity measure for Table 7
* local s "share_dummy"                                                                      // Alternative intensity measure for Table A3

qui foreach outcome in `outcomelist' {

use merged_base, clear

replace share=. if share==0

egen i1 = xtile(share) if year==2011,  n(4)
g i2 = (i1==4) if !missing(i1)
egen share_dummy = mean(i2) , by(fips)

keep if tin(2007,2016)
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

g treat_post = `s'*post

* With only year and county FEs (e.g., column 1 of Table 2)
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table7.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

* Add interaction between each county characteristic and year FE (e.g., columns 2 to 5 of Table 2)
foreach c in median_inc median_home_value education black_per {
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._`c' ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._`c') 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table7.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC
}

* Account for all controls (e.g., column 6 of Table 2)
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._median_inc i.year#c._median_home_value i.year#c._education i.year#c._black_per) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._median_inc i.year#c._median_home_value i.year#c._education i.year#c._black_per) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table7.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}


*** 2.1.3 Table A3 Alternative Treatment Intensity Specification

local s "share_dummy"                                                                      // Alternative intensity measure for Table A3

qui foreach outcome in `outcomelist' {

use merged_base, clear

replace share=. if share==0

egen i1 = xtile(share) if year==2011,  n(4)
g i2 = (i1==4) if !missing(i1)
egen share_dummy = mean(i2) , by(fips)

keep if tin(2007,2016)
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

g treat_post = `s'*post

* With only year and county FEs (e.g., column 1 of Table 2)
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using tablea3.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

* Add interaction between each county characteristic and year FE (e.g., columns 2 to 5 of Table 2)
foreach c in median_inc median_home_value education black_per {
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._`c' ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._`c') 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using tablea3.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC
}

* Account for all controls (e.g., column 6 of Table 2)
reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._median_inc i.year#c._median_home_value i.year#c._education i.year#c._black_per) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._median_inc i.year#c._median_home_value i.year#c._education i.year#c._black_per) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using tablea3.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}


*** 2.2 Dynamic regressions 

*** 2.2.1 Table 3 Dynamic Specification 

local s "treat_oil_gas"
local outcomelist " emp_all   avg_pay_all "

qui foreach outcome in `outcomelist' {

use merged_base, clear	

replace share=. if share==0

label var emp_all "EMP"
label var avg_pay_all "PAY"

keep if !mi(`s')
keep if tin(2007,2016)
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

forvalues i=2007/2016 {
	g yr_`i' = (year == `i')
	g treat_`i' = `s'*yr_`i'
	drop yr_`i'
	}

* With only year and county FEs (e.g., column 1 of Table 3)
reghdfe `outcome'_r treat_2007-treat_2016,  vce(cluster state)  absorb(i.fips i.year ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_2007-treat_2016, generate(n) absorb(i.fips i.year  ) 
forvalues i=2007/2016 {
	label var ntreat_`i' "`i'"
}
ols_spatial_HAC n`outcome'_r ntreat_2007-ntreat_2010 ntreat_2012-ntreat_2016, ///
lat(lat) lon(lng) t(year) p(fips) dist(200) bartlett lag(5) 
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'	
outreg2 using table3.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  label
drop n`outcome'_r-_est_spatHAC


* Add interaction between each county characteristic and year FE (e.g., columns 2 to 5 of Table 3)
foreach c in median_inc median_home_value education black_per {
reghdfe `outcome'_r treat_2007-treat_2016,  vce(cluster state)  absorb(i.fips i.year i.year#c._`c' ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_2007-treat_2016, generate(n) absorb(i.fips i.year i.year#c._`c' ) 
forvalues i=2007/2016 {
	label var ntreat_`i' "`i'"
}
ols_spatial_HAC n`outcome'_r ntreat_2007-ntreat_2010 ntreat_2012-ntreat_2016, ///
lat(lat) lon(lng) t(year) p(fips) dist(200) bartlett lag(5) 
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'	
outreg2 using table3.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  label
drop n`outcome'_r-_est_spatHAC

}

* Account for all controls (e.g., column 6 of Table 3)
reghdfe `outcome'_r treat_2007-treat_2016,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_2007-treat_2016, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education ) 
forvalues i=2007/2016 {
	label var ntreat_`i' "`i'"
}
ols_spatial_HAC n`outcome'_r ntreat_2007-ntreat_2010 ntreat_2012-ntreat_2016, ///
lat(lat) lon(lng) t(year) p(fips) dist(200) bartlett lag(5) 
matrix b = r(table)
matrix C_`outcome'_`s' = J(3,10,.)
matrix rownames C_`outcome'_`s' = coef ll95 ul95
matrix colnames C_`outcome'_`s' = 2007 2008 2009 2010 2011 2012 2013 2014 2015 2016
forvalues i = 1/4 {
matrix C_`outcome'_`s'[1,`i'] = b[1,`i'] \ b[5,`i'] \ b[6,`i']
}
matrix C_`outcome'_`s'[1,5] = 0 \ 0 \ 0
forvalues i = 6/10 {
local j = `i' - 1
matrix C_`outcome'_`s'[1,`i'] = b[1,`j'] \ b[5,`j'] \ b[6,`j']
}

sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'	
	
outreg2 using table3.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  label
drop n`outcome'_r-_est_spatHAC

}


*** 2.2.2 Figure 5 Treatment Effects of Coal Transition
coefplot matrix(C_emp_all_treat_oil_gas[1]), ci((C_emp_all_treat_oil_gas[2] C_emp_all_treat_oil_gas[3])) vertical title("Panel A: Employment") yline(0)  saving(f51, replace)  scheme(s2mono) ylabel(,angle(horizontal))  ylabel(, format(%9.1fc)) legend(off) 
coefplot matrix(C_avg_pay_all_treat_oil_gas[1]), ci((C_avg_pay_all_treat_oil_gas[2] C_avg_pay_all_treat_oil_gas[3])) vertical title("Panel B: Wages") yline(0)  saving(f52, replace)  scheme(s2mono) ylabel(,angle(horizontal))  ylabel(, format(%9.1fc)) legend(off) 
graph combine f51.gph f52.gph , scheme(sj) ysize(2) xsize(4) 
graph export fig5.tif, replace


*** 2.2.3 Figure A1 Treatment Effects of Coal Transition: 2002-2016

local s "treat_oil_gas"
local outcomelist " emp_all   avg_pay_all "

qui foreach outcome in `outcomelist' {

use merged_base, clear	

foreach v in   emp_all   avg_pay_all  {
 drop `v'_r
 g `v'_2002 = `v' if year==2002
 bys fips (year): carryforward `v'_2002, replace
 g `v'_r = 100*`v'/`v'_2002
 drop `v'_2002
}

keep if !mi(`s')
keep if tin(2002,2016)
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==15

forvalues i=2002/2016 {
	g yr_`i' = (year == `i')
	g treat_`i' = `s'*yr_`i'
	drop yr_`i'
	}

hdfe `outcome'_r treat_2002-treat_2016, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education ) 
forvalues i=2002/2016 {
	label var ntreat_`i' "`i'"
}

ols_spatial_HAC n`outcome'_r ntreat_2002-ntreat_2010 ntreat_2012-ntreat_2016, ///
lat(lat) lon(lng) t(year) p(fips) dist(200) bartlett lag(5) 

matrix b = r(table)
matrix C_`outcome'_`s' = J(3,15,.)
matrix rownames C_`outcome'_`s' = coef ll95 ul95
matrix colnames C_`outcome'_`s' = 2002 2003 2004 2005 2006 2007 2008 2009 2010 2011 2012 2013 2014 2015 2016
forvalues i = 1/9 {
matrix C_`outcome'_`s'[1,`i'] = b[1,`i'] \ b[5,`i'] \ b[6,`i']
}
matrix C_`outcome'_`s'[1,10] = 0 \ 0 \ 0
forvalues i = 11/15 {
local j = `i' - 1
matrix C_`outcome'_`s'[1,`i'] = b[1,`j'] \ b[5,`j'] \ b[6,`j']
}

}

coefplot matrix(C_emp_all_treat_oil_gas[1]), ci((C_emp_all_treat_oil_gas[2] C_emp_all_treat_oil_gas[3])) vertical title("Panel A: Employment") yline(0)  saving(fa11, replace)  scheme(s2mono) ylabel(,angle(horizontal))  ylabel(, format(%9.1fc)) legend(off) xlabel(,angle(45))
coefplot matrix(C_avg_pay_all_treat_oil_gas[1]), ci((C_avg_pay_all_treat_oil_gas[2] C_avg_pay_all_treat_oil_gas[3])) vertical title("Panel B: Wages") yline(0)  saving(fa12, replace)  scheme(s2mono) ylabel(,angle(horizontal))  ylabel(, format(%9.1fc)) legend(off)  xlabel(,angle(45))
graph combine fa11.gph fa12.gph ,  scheme(sj) ysize(2) xsize(4)
graph export figa1.tif,replace


*** 2.3 Alternative Control Groups

*** 2.3.1 Table 4 Alternative Control Groups

local outcomelist " emp_all  avg_pay_all "

use merged_base, clear	
rename (treat treat_100 treat_200 treat_psm) (Coal_States Within_100M Within_200M PSM)
tempfile temp
save `temp'

qui foreach outcome in `outcomelist' {
	
foreach s  in Coal_States Within_100M Within_200M PSM {

use `temp', clear

keep if tin(2007,2016)
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

g treat_post = `s'*post

reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & `s'==1 & year==2007
local dmean=`r(mean)'
outreg2 using table4.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}

}


*** 2.3.2 Table A2 Alternative Control Groups: Dynamic Specification

local outcomelist " emp_all  avg_pay_all "

use merged_base, clear	
rename (treat treat_100 treat_200 treat_psm) (Coal_States Within_100M Within_200M PSM)
tempfile temp
save `temp'

foreach s  in Coal_States Within_100M Within_200M PSM {

qui foreach outcome in `outcomelist' {

use `temp', clear	

keep if !mi(`s')
keep if tin(2007,2016)
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

forvalues i=2007/2016 {
	g yr_`i' = (year == `i')
	g treat_`i' = `s'*yr_`i'
	drop yr_`i'
	}
reghdfe `outcome'_r treat_2007-treat_2016,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education ) 
local r2 = 	`e(r2_a)'

hdfe `outcome'_r treat_2007-treat_2016, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education ) 
forvalues i=2007/2016 {
	label var ntreat_`i' "`i'"
}

ols_spatial_HAC n`outcome'_r ntreat_2007-ntreat_2010 ntreat_2012-ntreat_2016, ///
lat(lat) lon(lng) t(year) p(fips) dist(200) bartlett lag(5) 
sum `outcome' if e(sample) & `s'==1 & year==2007
local dmean=`r(mean)'	
	
outreg2 using tablea2.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  label
drop n`outcome'_r-_est_spatHAC

}

}


*** 2.4 Alternative outcome variables 


*** 2.4.1 Generate Table 5 

local s  "treat_oil_gas"
local outcomelist "emp_non_coal avg_pay_non_coal emp_fed avg_pay_fed emp_s_l avg_pay_s_l "  // For Table 5 Employment and Wages in Non-Mining and Public Sectors

qui foreach outcome in `outcomelist' {

use merged_base, clear

keep if tin(2007,2016)                            // For employment, wages, and population
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10                                // For employment, wages, and population

g treat_post = `s'*post

reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & `s'==1 & year==2007
local dmean=`r(mean)'
outreg2 using table5.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}


*** 2.4.2 Generate Panel A of Table 10

local s  "treat_oil_gas"
local outcomelist "popestimate pop_young pop_old pop_white pop_black "                      // Un-comment to generate Panel A of Table 10 

qui foreach outcome in `outcomelist' {

use merged_base, clear

keep if tin(2007,2016)                            // For employment, wages, and population
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10                                // For employment, wages, and population

g treat_post = `s'*post

reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & `s'==1 & year==2007
local dmean=`r(mean)'
outreg2 using table10.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}


*** 2.4.3 Generate Panel A of Table 11

local outcomelist "application_amt_1 application_amt_white_1 application_amt_black_1"       // Un-comment to generate Panel A of Table 11
local s  "treat_oil_gas"

qui foreach outcome in `outcomelist' {

use merged_base, clear

keep if tin(2008,2016)                             // For mortgage outcomes
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==9                                 // For mortgage outcomes 

g treat_post = `s'*post

reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & `s'==1 & year==2008
local dmean=`r(mean)'
outreg2 using table11.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}


*** 2.5 Table 6 Placebo Tests

forvalues i = 2000(3)2006 { 

use merged_base, clear	

foreach v in   emp_all   avg_pay_all  {
 drop `v'_r
 g `v'_`i' = `v' if year==`i'
 bys fips (year): carryforward `v'_`i', replace
 g `v'_r = 100*`v'/`v'_`i'
 drop `v'_`i'
}
drop post
local j = `i' + 5
local ii = `i' + 2
g post = (year>`ii')

tempfile temp
save `temp', replace

local outcomelist " emp_all  avg_pay_all  "
local s "treat_oil_gas"

qui foreach outcome in `outcomelist' {

use `temp', clear	

keep if tin(`i',`j')
keep if !mi(`s')
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==6

g treat_post = `s'*post

reghdfe `outcome'_r treat_post ,  vce(cluster state)  absorb(i.fips i.year ) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat_post, generate(n) absorb(i.fips i.year ) 
ols_spatial_HAC n`outcome'_r ntreat_post , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) lag(5)
sum `outcome' if e(sample) & `s'==1 & year==`i'
local dmean=`r(mean)'
outreg2 using table6.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'_r-_est_spatHAC

}

}


*** 2.6 Access to finance and mobility

*** 2.6.1 Generates Table 8 and Table 9, and Panel B of Table 10

local outcomelist "emp_all avg_pay_all popestimate"     

local i = 8

qui foreach outcome in `outcomelist' {

use merged_base, clear

replace treat = treat_oil_gas

keep if !mi(treat)
keep if tin(2007,2016)
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==10

rename (low_atf_bra low_large_bra low_mob_rel)  (banks large  mob)

* Columns (1) to (4)
foreach a2f in banks large rural mob {
reghdfe `outcome'_r treat#post post#`a2f'  treat#post#`a2f',  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat#post post#`a2f'  treat#post#`a2f', generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r n__1__* , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) 
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table`i'.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'-_est_spatHAC

}

* Columns (5) to (7)
foreach a2f in banks large rural  {
reghdfe `outcome'_r treat#post  post#`a2f'  treat#post#`a2f'   ///
 post#mob  treat#post#mob,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
hdfe `outcome'_r treat#post  post#`a2f'  treat#post#`a2f'   ///
 post#mob  treat#post#mob, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r n__1__* , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) 
sum `outcome' if e(sample) & treat==1 & year==2007
local dmean=`r(mean)'
outreg2 using table`i'.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'-_est_spatHAC
}

local i = `i' + 1

}


*** 2.6.2 Generate Panel B of Table 11

local outcomelist "application_amt_1  "    

qui foreach outcome in `outcomelist' {

use merged_base, clear

replace treat = treat_oil_gas

keep if !mi(treat)
keep if tin(2008,2016)
drop if missing(lat,lng,`outcome'_r)

gen temp=1 if `outcome'_r!=. 
bys fips: egen num_obs=total(temp)
keep if num_obs==9

rename (low_atf_bra low_large_bra low_mob_rel)  (banks large  mob)

* Columns (1) to (4)
foreach a2f in banks large rural mob {
reghdfe `outcome'_r treat#post post#`a2f'  treat#post#`a2f',  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
local r2 = 	`e(r2_a)'
hdfe `outcome'_r treat#post post#`a2f'  treat#post#`a2f', generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r n__1__* , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) 
sum `outcome' if e(sample) & treat==1 & year==2008
local dmean=`r(mean)'
outreg2 using table11.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'-_est_spatHAC

}

* Columns (5) to (7)
foreach a2f in banks large rural  {
reghdfe `outcome'_r treat#post  post#`a2f'  treat#post#`a2f'   ///
 post#mob  treat#post#mob,  vce(cluster state)  absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
hdfe `outcome'_r treat#post  post#`a2f'  treat#post#`a2f'   ///
 post#mob  treat#post#mob, generate(n) absorb(i.fips i.year i.year#c._black_per  i.year#c._median_inc i.year#c._median_home_value  i.year#c._education) 
ols_spatial_HAC n`outcome'_r n__1__* , ///
lat(lat) lon(lng) t(year) p(fips) dist(200) 
sum `outcome' if e(sample) & treat==1 & year==2008
local dmean=`r(mean)'
outreg2 using table11.doc,  append $options addstat(R2, `r2', Dependent Variable Mean, `dmean')  
drop n`outcome'-_est_spatHAC
}

}


forvalues i = 2/11 {
	erase table`i'.txt
}

foreach v in a2 a3 {
	erase table`v'.txt
}

foreach v in f31 f32 f33 f34 f41 f42 f43 f44 f2 f51 f52 fa11 fa12 {
	erase `v'.gph
}
